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ABSTRACT 



We report the discovery of kHz fluctuations, including quasi-periodic oscillations (QPO) 
at ~ 330 Hz (260 — 407 Hz) and ~ 760 Hz (671 — 849 Hz) and a broadband kHz continuum in 
the power density spectrum of the high mass X-ray binary pulsar Centaurus X-3 ( |Jernigan^ 
Klein and Arons 1999| ). These observations of Cen X-3 were carried out with the Rossi 
X-ray Timing Explorer (RXTE). The fluctuation spectrum is flat from mHz to a few Hz, 
then steepens to f~ 2 behavior between a few Hz and ~ 100 Hz. Above a hundred Hz, 
the spectrum shows the QPO features, plus a flat continuum extending to ~ 1200 Hz and 
then falling out to ~ 1800 Hz. These results, which required the co-adding three days of 
observations of Cen X-3, are at least as fast as the fastest known variations in X-ray emission 
from an accreting compact object (kHz QPO in LMXB sources) and probably faster since 
extension to ~ 1800 Hz is indicated by the most likely parameterization of the data. 

Multi-dimensional radiation hydrodynamics simulations of optically thick plasma flow 
onto the magnetic poles of an accreting neutron star show that the fluctuations at frequencies 
above 100 Hz are consistent with photon bubble turbulence and oscillations (PBO) previously 
predicted ([Klein et al. 1996a|) to be observable in this source. We show that previous 
observations of Cen X-3 constrain the models to depend on only one parameter, the size of 
the polar cap. For a polar cap opening angle of 0.25 radians (polar cap radius ~ 2.5 km and 
area ~ 20 km 2 , for a neutron star radius of 10 km), we show that the spectral form above 100 
Hz is reproduced by the simulations, including the frequencies of the QPO and the relative 
power in the QPO and the kHz continuum. This has resulted in the first model-dependent 
measurement of the polar cap size of an X-ray pulsar. The simulations underpredict the 
overall amplitude of the observed spectrum, which we suggest is the consequence of a 2D 
axisymmetric simulation of an intrinsically 3D phenomenon. The power density spectrum of 
Cen X-3 shows a dramatic decrease above ~ 1000 Hz which suggests an optical depth ~ 30 
across the accretion mound consistent with effects of radiative diffusion in the simulations. 
We identify this decline at high frequency as the first direct evidence of radiative diffusion 
near the surface of a neutron star (NS). 

We suggest the fluctuations observed at frequencies below 100 Hz, whose spectrum has 
a different form from that of the kHz phenomena, reflect intermittency in the mass transfer 
mechanism which carries plasma from the accretion disk to field aligned flow onto the neutron 
star's polar caps. Using simple estimates based on Rayleigh- Taylor instabilities, possibly 
modulated by intrinsic disk turbulence, we show that mass transfer in "blobs" forming 
through Rayleigh- Taylor disruption of the disk's inner edge can explain the large amplitude 
fluctuations required by the spectrum at frequencies / ~ 1 Hz, but only if magnetic pressure 
in the disk's innermost regions inflates the disk until its scale height is comparable to the 
magnetosp here's size ~ 4300 km. 

The observational results required the development of a procedure for the careful de- 
termination of the deadtime effects of the PCA. This procedure is described in appendix 
A. As a consequence of the use of observations of Cyg X-l for the estimation of deadtime 



2 



corrections for the observations of Cen X-3, we have also demonstrated that the black hole 
Cyg X-l shows clear evidence of variability up to a frequency of ~ 280 Hz. Also observations 
of GX 17+2 were used to validate the procedure for deadtime corrections. This analysis of 
GX 17+2 clearly indicates the presence of a kHz QPO and the absence of any significant 
simultaneous kHz continuum. 



1. Introduction 



The binary system composed of the X-ray emitting neutron star Cen X-3 and its companion O-type 
supergiant has been studied extensively by all orbiting X-ray astronomy observatories. The results from 
the Uhuru satellite identified the system as an eclipsing binary X-ray pulsar with a 2.087 day orbit and 
a 4.84 s pulse period ( Giacconi et al. 1973 ; Schreier et al. 19?"2|) . The details of the orbit determination, 
which include estimates of the masses of both the neutron star and its companion, are summarized in 
a review ( [Rappaport and Joss 198^ ). The well determined distance to the binary system ( [Hutchings ei 
al. 1979| ) yields an accurate determination of the luminosity. Recently a cyclotron feature was detected 



in the X-ray spectrum of Cen X-3 providing precise knowledge of its surface magnetic field strength 
( [Santangelo at al. 1998| ). 
Table 1 summarizes the parameters of Cen X-3. The average X-ray luminosity is 



L x = 9.4 x 10 37 



D 



10 kpc 



Z) 



1.3 



erg s 



where Z = ((1 — 2GM/Rc 2 ) x l 2 — 1) is the surface gravitational redshift of the star of mass M and 
radius R. The surface magnetic field strength B implied by the observed cyclotron feature is 



B 



3.2 x 10 



12 



E c 



28.5 keV 



1.3 



Gauss 



(2) 



where E c is the energy of the cyclotron line. The measured mass of the neutron star combined with 
theory for the internal structure ( Baym and Pet hick 1979| ) constrains the range for its radius. If one 
assumes the magnetic field geometry at the surface to be dipolar, then the size of the polar cap onto 
which matter accretes is not well constrained by direct measurement. 
EDITOR: PLACE TABLE 1 HERE. 

Motivated by earlier predictions of the likelihood of the existence of Photon Bubble Oscillations (PBO) 
in Cen X-3 ([Klein et al. 1996a| ), we have analyzed two consecutive binary cycles of the source to search for 
PBO. We have also calculated new theoretical models constrained by the observable parameters of Cen 
X-3 (Table 1). The value for L x is consistent with the observations reported here. The physics describing 
the radiation hydrodynamics which governs the accretion of matter onto the highly magnetized polar 
caps of luminous X-ray pulsars, and which has been incorporated in these new numerical simulations, 
has been described elsewhere ( |Arons, Klein, and Lea 1987| ; [Klein and Arons 1989| ). Early numerical 
results ( [Klein and Arons 1991| ) and linear stability analysis ( |Arons 1992| ) suggested the formation of 
small scale but large amplitude fluctuations in the matter density and velocity, which form radiation 
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filled pockets almost devoid of plasma ("photon bubbles"). These photon bubbles (PB) are embedded in 
the optically thick, inflowing plasma and as they grow non-linearly, they result in significant observable 
fluctuations in the emitted luminosity. 

We observed the X-ray emission from Cen X-3 with the Rossi X-ray Timing Explorer (RXTE) over 
two binary cycles of the source and clearly detect both kHz QPO and a broadband kHz continuum in 
the power density spectrum which we identify with PBO and the high frequency power law continuum 
generated by photon bubble "turbulence" . (Our use of the descriptive noun turbulence should not be 
confused with the classical turbulence of a fluid.) Given the observationally derived physical parameters 
for Cen X-3 which are summarized in Table 1, the only remaining free parameter in our numerical models 
is the size of the polar cap. In this paper, we show that by adopting a reasonable size for the polar cap 
of Cen X-3, we are able to semi- quantitatively match the frequencies of two observed QPO/PBO peaks 
and the kHz continuum with the power density spectrum of our calculated radiation-hydrodynamic 
models; thus giving strong evidence for the existence of PBO and photon bubble turbulence in Cen X-3. 



2. Observations 

Separate continuous observations of Cen X-3 during 1997 February-March (JD2450507.5067-JD2450509.2248; 
JD2450509.5470-JD245051 1.3304) were carried out with the PCA on-board RXTE ( Bradt, Rothschild, 



and Swank 1993| and |Swank 1997| ) for two consecutive binary cycles of the source. Each observation 



lasted 1.7 days; during the non-eclipse phases of the binary orbit, data were acquired in two broad-band 
energy channels at a time resolution of 62 /zs. The observation modes were specially selected to provide 
detailed information on the high time resolution behavior of the 4.8 s pulse as well as any broadband 
temporal behavior at frequencies above ~ 100 Hz. Figure |TJ is a log-log presentation of the average 
deadtime corrected power density spectrum for the entire data train of the low energy channel (1-7 
keV). The curve was computed by co-adding ~ 3400 separate power density spectra, each with a ~ 65 
s exposure and each composed of 2 20 bins of duration 62 fis. These power density spectra were further 
averaged over frequency bins evenly spaced in logarithmic steps. This approach is the most sensitive 
way to search for broadband features in the power density spectrum over a wide range of frequencies 
(10~ 2 — 10 4 Hz). Note the smooth plateau at very low frequencies (10~ 2 — 10 _1 Hz) followed by the first 
and second harmonics of the 4.8 s pulse (10 _1 — 10° Hz). The significance of the pulse is degraded by 
the crude averaging over frequency. However, the integral of the power density still correctly indicates 
the large pulsed fraction of the source. The periodic pulse and its harmonics are superimposed on a 
broadband power law which falls from an amplitude of ~ 10 _1 to 10~ 6 with a slope of approximately —2 
(see dashed line). At about 100 Hz another power density plateau appears, followed by two QPO peaks 
(~ 331 and ~ 761 Hz). These two QPO peaks (50% confidence ranges 260 - 407 Hz and 671 - 849 Hz) 
are superimposed on the power density plateau which falls sharply at frequencies above ~ 1200 Hz. The 
curve falls to a level of about 10 -7 ( -2S§_] Hz -1 at ~ 1800 Hz which has reached the Poisson detection 

\ mean J 

limit of the transform. The high frequency Poisson limit in the form a frequency dependent function 
has been subtracted from the transform of the data (equation 6 of Appendix A). The Leahy normalized 
version of the transform ( Leahy et al. 198~3|) was carefully corrected for the effects of deadtime of the 
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PCA (see Appendix A). 

The same analysis procedure was applied to the high energy channel (> 7 keV). For these data the 
power density spectrum is nearly identical to that shown in Figure [I] out to a frequency of about 50 
Hz. We detect no significant kHz features. Determination of the energy dependence of the kHz features 
of the power density spectrum must await significant additional PCA data or future detections with 
instruments of much larger aperture than RXTE. 

EDITOR: PLACE FIGURE [I] HERE. 

We note that the identification of the specific features in the power density spectrum as QPO in the 
presence of a significant continuum is a somewhat arbitrary description. Others may prefer to identify 
the QPO as "bumps" in the continuum or just as a complex continuum. In this paper, we choose to use 
the QPO designation in analogy with its use for the discovery of QPO in the LMXB source GX5-1 flyan 
der K jis et ai. 1995|Jvan der Kris et al. 1995] ). We justify this view based on the similar amplitude and 
width of these QPO in the HMXB Cen X-3 and the low frequency (~ 6 Hz) QPO in the LMXB GX5-1. 
Previous work on QPO in X-ray sources has indicated a large dynamic range of Q values including 
examples of order unity. The formal significance of the detection of these kHz QPO in Cen X-3 is not 
as high as for the QPO in GX5-1 but this fact is likely just a coincidental result of the limited aperture 
of the PCA on-board RXTE. Even if these kHz features are real statistical fluctuations in behavior of 
Cen X-3, the QPO designation is still an appropriate description of the results presented in this paper. 
In using the QPO designation by observational analogy with GX 5-1 we do not suggest any underlying 
physical analogy. 

A model composed of fifteen free parameters was fitted to the power density over the frequency range 
40 — 4000 Hz. These parameters are a low corner frequency Fi, a high corner frequency F^, a level of 
power density P[ at the frequency i 7 } and its constant extension to lower frequencies, the slope Si of 
a power law below the frequency F[ and its extension to higher frequencies, the slope S m of a power 
law between the frequencies Fi and Fh, the slope Sh of a power law above the frequency Fh and its 
extension to higher frequencies, the frequency of the first Gaussian QPO Fx, its rms amplitude A\ and 
its width Wx-i and the frequency of the second Gaussian QPO F 2 , its rms amplitude A 2 and its width 
W%. Finally the last three parameters (A, B and C) model the frequency dependent deadtime correction 
to the power density spectrum (see Appendix A for details). 

Alternative forms of the continuum power density will slightly change the values of the six parameters 
that describe the two QPO peaks but will not change the basic interpretation of the fit. This particular 
functional form for this model was selected as a reasonable choice for a simple characterization of the 
full power density spectrum between 40 — 4000 Hz with particular attention to the frequency range 
~ 100 - 2000 Hz. 

Table 2 displays a list of the parameters and their best fitted values, as determined by a minimization 
of x 2 - The errors in the measured power density points were determined by the Poisson counting noise. 
These data are comprised of 39 independent points, resulting in a minimum x 2 per degree of freedom of 
1.38 for 24 degrees of freedom. Formally this is not an adequate fit, however inspection of the residuals 
indicates that the fit is within expected limits for frequencies between 100 — 2000 Hz (see Figure |9| in 
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Appendix A). 

All components of the model match the data well except for small deviations at low frequencies (see 
Figure |2| and Figure ^|). The minor deviations from a power law at low frequencies (< 100 Hz) are 
revealed by the extremely small error bars. These deviations which are a small fraction of the power 
density at low frequency are likely caused by a real physical mechanism that is more complex than a 
simple power law. Since the primary goal of this paper is to characterize the kHz behavior of Cen X-3 
we did not further complicate the fifteen parameter model. The high frequency data are of insufficient 
signal to noise to realistically determine the exact nature of this detailed structure in the falling power 
density. The uncertainty in the interpretation of the decline in power above ~ 1200 Hz is quantified by 
the large acceptable range of the parameter Sh- The presence of the two QPO peaks and the broadband 
kHz continuum in the power spectrum over the frequency range 100 — 1800 Hz is highly significant. 
Settting the amplitudes of the QPO peaks to zero results in an increase in x 2 of ~ 47 and ~ 42 for QPO 
1 and QPO 2 respectively. This implies that the first and second QPO peaks are formally detected 
at greater than 99.9% confidence. The lower and upper values of the parameters in table 2 are 50% 
confidence limits for each of the fifteen free parameters of the model. (These are determined by the 
variation of each parameter until y 2 increased by ~ 14.34, see [Lampton, Margon and Bowyer 1976| .) 



The form of the fitted model for the low frequency (40 — 100 Hz) power law was selected with the 
assumption that whatever the physical mechanism that forms this power law that behavior is likely 
to continue past 100 Hz. Since the extrapolation of the low frequency power law above 100 Hz is 
much lower than the observed power density in the kHz band this suggests a transition near 100 Hz in 
which a new physical mechanism becomes dominant. Note that the continuation of the power law to 
very low frequency (< 3 Hz) rolls over to a flat level as required by causality since this extrapolation 
cannot produce fluctuations that exceed 100% rms. This damping of the behavior at low frequency is 
still consistent with a single physical mechanism. Possibilities for the physics that generates this low 
frequency power law are discussed in a later section. 

EDITOR: PLACE TABLE 2 HERE. 

Figure ^| shows the deadtime corrected data with ±1<t error bars (dotted lines) as well as the model 
(solid line). The data have been multiplied by frequency in order to more clearly indicate the Poisson 
counting noise. This also means that the vertical axis has units of (rms/mean) 2 which allows for the 
easy approximate interpretation of features in units of rms percentage amplitude. Dashed horizontal 
lines at 2, 5 and 100% are shown to indicate the relative amplitudes of features of the power density 
spectrum. The coherent pulsations from Cen X-3 clearly have a large amplitude of order ~ 50% whereas 
the high frequency components (> 100 Hz) have amplitudes in the 2-5% range. The best fitted fifteen 
parameter model is also plotted in Figure § and clearly tracks the measured data. Note the presence of 
the two kHz QPO peaks (~ 330 Hz and ~ 760 Hz) and the break in the continuum at the high corner 
frequency (F h ~ 1200 Hz). 

EDITOR: PLACE FIGURE || HERE. 
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3. Theoretical Model 

Photon bubbles instabilities can form in magnetized radiation pressure supported atmospheres and 
subsonic flows ( [Klein and Arons 1989| , |Arons 1992| ) with B fields as small as 10 8 Gauss. Classical X-ray 



pulsars such as Cen X-3 easily fulfill this criterion. We calculate a model appropriate for the surface 
conditions of the neutron star Cen X-3. This model relies on the self-consistent solution of the full two- 
dimensional, time-dependent equations of the radiation hydrodynamics governing the accretion of matter 
onto the highly magnetized polar caps of a luminous x-ray pulsar QArons, Klein, and Lea 1987 , [Klein and 



Aron sp~989| , |Klein et al. 1996a| ). The full description of the physics of highly magnetized flow that we 



have included in our calculations is described in detail in Arons, Klein, and Lea 1987 . We solve the full 
set of two-dimensional, time-dependent relativistic radiation hydrodynamic equations. These include 
the equations which govern the conservation of mass and momentum, and separate energy equations for 
both the electrons and the ions as well as a full non-LTE description of the time dependent transport of 
radiation and photon number density ( [Klein and Arons 1991| , and [Klein et al. 1996c| ). We also include 
the full effects of super strong magnetic fields on all physical processes and we have developed a new 
flux-limited diffusion theory for the transfer of radiation in optically thin, radiation filled structures 
embedded in optically thick flows (Photon Bubbles). We solve the entire structure of self-consistent 2D 
coupled non-linear P.D.E's with a 2nd-order accurate, time explicit operator splitting approach. 

The full evolution of the radiation hydrodynamic equations is an initial value problem in time. The 
surface value of the magnetic field strength assumed to have dipole geometry, the mass accretion rate, 
and the size of the polar cap are all parameters, as are the mass and radius (which enter into the strength 
of gravity) of the neutron star. We hold these characteristics fixed in the course of a calculation (see 
Table 1 for values). This leaves the geometry of the polar cap as the only free parameter in the 
calculation. In principle, a specification of the surface field strength B and the intrinsic luminosity L x 
could, under certain assumptions (i.e., magnetic dipole field configuration) determine the radius of the 
polar cap 6 C ( |Pringle and Kees 1972| , |Davidson and Ustriker 1973| , [Lamb, Pethick and Pines 1973] , [Arons 



and L ea 198(f ). The details of the magnetospheric mass loading are quite uncertain. Since these details 



determine the size of the polar cap and the distribution of the accretion mass flux within the polar flux 
tube, we simplify by assuming the accretion flux to be constant from the magnetic axis out to a fixed 
magnetic flux surface, then zero beyond, and we take the polar cap radius to be the only free parameter 
in our calculations. 

From our previous calculations with L x = 3 x 10 37 ergs s _1 and B = 3 x 10 12 G, we know that in models 
with large polar caps (9 C = 0.4) we obtain two PBOs with frequencies of 170 Hz and 260 Hz ( Klein ei 



al. 1 9,96b| ) while smaller polar caps (6 C = 0.1) yield PBO with significantly higher frequencies in the 



range of 3500 Hz ( |Klein et al. 1996a]) . Thus wider accretion columns produce PBO of lower frequency, a 
scaling most likely related to the longer diffusion time necessary for photons which originate in interior 
photon bubbles of wider columns to diffuse to the outer edge of the accretion column. Therefore, PBO 
at the observed frequencies at 330 Hz and 760 Hz should occur for polar cap sizes bounded by our 
previous calculations, 0.1 < 9 C < 0.4. We have computed three new models to explore the variation of 
PBO frequency with polar cap sizes 8 C = 0.25, 0.3 and 0.4, with the luminosity L x = 4.7 x 10 37 ergs 
s _1 and magnetic field strength B = 3.2 x 10 12 G inferred from observations of Cen X-3. The value of 
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L x selected here is one half of that shown in Table 1 which is appropriate for one of the two polar caps. 
The calculations refer to accretion onto a single polar cap. 

Since the basic oscillation and instability time scales of the magnetically modified g-modes which be- 
come photon bubbles are proportional to the time scale for radiation diffusion ( Arons 1992| ), and photon 



escape from the sides of the column controls the emergent intensity, one can reasonably suppose that 
the time scale for radiation diffusion across the magnetic field, — -^i(3/cA m / Pi j_) = L±/Vdiff± oc 9\, 
might be related to PBO timescales. This strong dependence on 9 C is amply supported by the simula- 
tions previously published and reported here. Here V<nff± = c/3tj_ is the diffusion velocity across the 
flow above the accretion shock, with t± the scattering optical depth across the magnetic field. L± is the 
transverse distance across the magnetic field which at the base of the accretion mound on the stellar 
surface is -R*^, but further up in the accretion column is i?*^ c (r/i?*) 15 . The parameter A m / Pi j_ is the 
transverse mean free path across the accretion mound. We find in the 9 C = 0.25 model that just above 
the accretion shock at the accretion mound, tj_ = 19 and L± = 6 x 10 5 cm, yielding — 1.14 x 10~ 3 s. 
This corresponds to a diffusion frequency Udiff — 877 Hz. At the base of the accretion column, we find 
t_i_ = 80 yielding ^// = 2 X 10~ 3 s and the corresponding diffusion frequency Vdiff = 500 Hz. 

Calculating the transverse radiation transport speed as F±/E rad at the accretion shock, where F± is the 
transverse radiation flux and E ra d is the radiation energy density in the 9 C = 0.25 model, and setting the 
result equal to Vdiff-, gives a diffusion speed of 1.0 x 10 9 cm s -1 and a diffusion time tdiff = 6.04 x 10~ 4 



'diffi gives a uiiiusiuu speeu ui i.u a iu cm s tuiu a uniusiuu iime tdiff 
s. This corresponds to a diffusion frequency Pdiff = 1-6 X 10 3 Hz. Taking the measurement of the 
transverse radiation transport speed near the base of the accretion column gives Vdiff — 5 x 10 8 cm 
s~\ tdiff = 5 x 10~ 4 s and Vdiff = 2000 Hz. Thus the two ways of computing the diffusion time are 
in agreement within roughly a factor of 2. Considering these approximations, one might expect this 
model to show PBO with fundamental frequency fpso ~ 1000 Hz, not drastically different from the 
values found in the computational model and in the observations of Cen X-3. The uncertainty in the 
quantification of the diffusion time to within a factor of 2 is a consequence of the highly non-uniform 
structure of the accretion column. If the column were uniformly dense, the two approximations for Tdiff 
would be in more precise agreement. This factor of 2 uncertainty has no bearing on the accuracy of the 
PBO centroid frequencies which are determined by the solutions of the full radiation hydrodynamics 
equations in a spatially inhomogenious medium. Thus the final PBO centroid frequencies represent 
a complex path for photons through intermingled optically thick and thin plasma. A diffusion time 
estimated at a few points in such a medium is not likely to be better than a factor of 2 in accuracy. 
A Monte Carlo treatment following the time history of the photons in the accretion column would be 
necessary to get a precise determination of the diffusion time for photons to leak out of the column. 

We identify the high corner frequency ~ 1200 Hz noted previously with an average diffusion time 
scale tdiff which we relate to the total optical depth r & s across the accretion mound. 



r obs = 26 




1 


9 C 


-l 


r r* - 


-l 




0.25 radians 




.10 km_ 





(3) 



This equation for T b. s is a recasting of the expression for the diffusion time tdiff, now expressed in terms 
of the optical depth. The form of Equation 13] indicates how the value of T b s varies depending on the 
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assumed size of the polar cap, 6 C . The particular value of 9 C is the assumed value that we used for model 
A which also corresponds to a typical value within the expected range for 8 C (see Table 1). Basically 
any internal oscillations within the interior of the accretion mound which are faster than the high corner 
frequency are suppressed by a factor that scales approximately as the inverse second power of frequency. 
This simple scaling is the result of the approximation that diffusion is just a smoothing process by a 
linear filter with a characteristic time scale. The estimated value of the optical depth r {, s ~ 26 agrees 
approximately with the values of optical depth in the numerical calculations (typically ~ 20 — 80). This 
agreement is not expected to be better than a factor of a few because of the simplifying assumption of 
a uniform medium with a diffusion velocity (3^-)- We identify this feature at high frequency as the 
first direct evidence of radiative diffusion near the surface of a neutron star (NS). 



4. Results of Radiation Hydrodynamic Calculations 

Figure shows the results for the power density spectra of the time series of the emergent transverse 
luminosity out of the sides of the accretion column for three models (A, B and C) for which the values 
L x and B have been set to those measured for Cen X-3 (see Table 1). The models A, B and C 
are distinguished by setting the polar cap size to 0.25, 0.3 and 0.4 radians respectively. The power 
density spectrum for each model is obtained by first performing a full 2-D time-dependent radiation- 
hydrodynamic calculation and obtaining the radiative flux in the transverse direction. This flux is 
integrated along a series of annular rings along the sides of the accretion funnel to obtain the emergent 
luminosity which exits the accretion column in a low altitude fan beam of radiation. The time series of 
the transverse luminosity is then Fourier transformed to obtain the power density spectrum. 

EDITOR: PLACE FIGURE | HERE. 

These models are simulations of only a few tens of milliseconds of accretion onto the neutron star 
and may not be well defined ergodic samples necessary for an exact comparison with the observed 
power density spectra averaged over many thousands of seconds. However, in support of making such 
comparisons of simulations to observations, we point out that any initial transients in the calculations 
die out after less than 10% of the running time for each simulation. The power spectra calculated from 
the simulations are derived from the remaining ~ 90% of each run, corresponding to what appears to 
be a statistically steady state. Our procedure limits the quality factor, Q, of any individual PBO to the 
ratio of the simulation's length to a given PBO period (Q < 100). Our procedure obviously would miss 
slow changes of the apparent statistically steady states observed. We have no reason to believe that 
there are such slow systematic effects, given the physics incorporated in the simulations. Therefore, 
we take the simulations as our best known estimates of the behavior of PBO over times much longer 
than characteristic PBO periods. We further note that the errors on the actual calculation of the power 
densities are very small (less than 10~ 8 ) in the units shown in figure |3]. The detailed features in the 
power densities are formally very significant with none of the degradation at high frequencies that typify 
observations containing Poisson counting noise. These very small errors, which cannot practically be 
shown in Figure [3|, are due to the finite spatial grid and discrete time steps of the numerical calculations. 
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This level of this numerical noise in the power density is revealed by an artificial downturn at very high 
frequencies (~ 10 5 Hz). Typically, timesteps in the calculations are < 5 x /vpbo- 

In Figure |3|, the lower curve is model A (8 C = 0.25 radians), the middle curve is model B (8 C = 0.3 
radians) and the upper curve is model C (8 C = 0.4 radians). Models B and C have been artificially shifted 
by two and four decades vertically for illustrative purposes. The curves in the left panel are shown in 
high resolution to reveal the details of structure of the power density spectrum. Below a frequency of 
~ 1 kHz the spectrum is shown at the full resolution of the calculations. At higher frequencies the data 
are smoothed in equal logarithmic intervals for clarity of presentation. The curves in the right panel 
show the same data further smoothed to facilitate comparison with observations (see figure 0) and to 
enhance the broadband features. The most striking common feature of all three models is the rolloff at 
high frequencies in approximate agreement with Equation |3| and a clear demonstration of the signature 
of radiative diffusion. In detail Models A and B show a somewhat steeper slope at higher frequencies 
which is not are that unlikely considering the rough approximation of Equation |3|. Each curve has a 
principal low frequency PBO that is shifted to higher frequency for the models with smaller polar caps in 
good agreement with vpbo 1/0^ scaling. The lower frequency PBO in model A determined from the 
detailed calculations is ~ 350 Hz (not redshifted), which agrees well with the observed lower frequency 
QPO for Cen X-3 (~ 330 Hz). Our detailed calculations yield lower frequency PBO approximately at 
250 Hz and 200 Hz for models B and C with 6 C = 0.30 radians and 9 C = 0.40 radians respectively. These 
are in reasonable agreement with the frequencies 246 Hz and 138 Hz suggested by d\ scaling. Also the 
level of the continuum is somewhat lower for models with smaller polar caps. 

Note that the absolute levels of the power density can be directly compared to the observed levels 
in Cen X-3 (see Figure [J). The amplitude of the observed power density in the range 100 — 1800 Hz 
is about one or two orders of magnitude higher than the level calculated in the theoretical models. 
This discrepancy is not that severe considering the simplified model for the geometry of the emission 
region. The actual rms amplitudes of the fluctuations in the models near 1 kHz are only a factor of ten 
lower than the observed level. The variability of the escaping X-ray luminosity is a damped form of the 
more violent internal photon bubble variability deep within the accretion mound. The presence of the 
PBs within the mound necessarily implies the local internal opacity variations have rms amplitudes that 
reach nearly 100%. In an extension of the 2D theoretical models to 3D we would expect the variations of 
opacities "seen" by internally generated photons to vary by an even wider range. Some of these photons 
might "see" reduced opacities episodically which would significantly enhance the X-ray variability and 
rms amplitudes seen by a distant observer. Also the assumption of a uniformly filled accretion column 
might be relaxed as a partially filled column (hollow cone) would also provide lower opacity paths from 
the interior of the shock mound to the exterior surroundings. We are extending our calculations to a 
hollow cone accretion geometry as the first step toward investigating 3D effects on the amplitude of the 
PBO fluctuations. This work will be reported elsewhere. There are likely very low frequency variations 
(< 1 s or pulsed modulated at 4.8 s) in the geometry of the matter filling the accretion column that 
would result in low opacity "windows" which would allow a distant observer to "see" into the interior 
of the shock mound where the rms fluctuations of the X-ray emission are much higher. 
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Model A (6 C = 0.25 radians) shows PBO at 350 Hz and a composite PBO at about 700 Hz. The 
frequencies of these PBO are in good agreement with frequencies of the QPO peaks observed in the 
power density spectrum of Cen X-3. The relative amplitudes of the PBO/QPO as compared to the 
continuum power are a better match between the theory and observation than the absolute amplitudes. 
The rms amplitudes of the observed QPO 1 and 2 in Table 2 can be expressed as equivalent widths 
(~ 710 Hz and ~ 970 Hz) which are comparable to the full range of the continuum (~ 100 Hz to 
~ 1800 Hz) indicating an approximately equal contribution of the QPOs and the continuum to the full 
rms amplitude above ~ 100 Hz. The approximate balance between PBO/QPOs and continuum is also 
roughly true for all three theoretical models. The additional strong broad composite PBO at ~ 3000 Hz 
that is seen in models A and B is not seen in the observed power density spectrum of Cen X-3. However 
it is below the level that could be easily detected with the PCA. The upper limit for the detection of 
broad features in the power density spectrum rises proportional to frequency above ~ 1000 Hz which 
rapidly decreases the capability for detection of any features near 3000 Hz. This rapidly declining ability 
to detect broad features at high frequencies should not be confused with the flat asymptotic behavior 
for the detection of features of fixed bandwidth that have increasing Q values at higher frequencies. 

Models A, B and C have been plotted simultaneously with the observed data in Figure |4]. For 
comparison both models and the observed data have been plotted over an appropriate frequency range 
and on the same absolute scale with no effort to adjust the amplitudes to match the data. The frequencies 
of the theoretical simulations have been red shifted correctly by the required factor of (1 + Z) = 1.3 
appropriate for the surface of a neutron star with mass and radius specified in Table 1. The frequencies 
of the PBO in model A match the frequencies of the QPO in the observed data. This is remarkable 
considering that the only free parameter in our calculations is the size of the polar cap. The theoretically 
determined PBO at ~ 300 Hz and the corresponding observed QPO at a similar frequency are both 
relatively weak in comparison to the PBO/QPO at ~ 700 Hz. In addition model A (also B) has a 
broad PBO with a large amplitude at ~ 3000 Hz that is very near the limit of the observable range of 
RXTE (Figure |3|). The continuum shape of model A does not match the observations as well as Model 
C (Figure |J). All three models show a power law rolloff at high frequency (~ 10 4 Hz) in agreement with 
the simple approximation that diffusion is just a smoothing process by a linear filter with behavior f~ 2 . 
Models A and B show a somewhat faster rolloff at even higher frequencies. If we consider the whole 
picture, model A has the best QPO/PBO frequency match and models B or C have a better match to 
the continuum. These models taken together determine a range for the polar cap radius of ~2-3 km. 
Calculations are underway (to be reported elsewhere) to better bracket the polar cap radius by doing 
a parameter study of small variation of the polar cap radius around our best model A with 9 C = 0.25 
radians, and its effect on the PBO frequencies. These calculations are extremely time consuming and 
typically require ~ 300 hours per model on a supercomputer. 

EDITOR: PLACE FIGURE | HERE. 

Our calculations have used measured parameters for Cen X-3 (e.g., luminosity, magnetic field, mass, 
etc.) while only varying the size of the polar cap. Hence the identification of the kHz QPOs and a 
kHz broadband continuum with PB phenomena suggests that the radius of the polar cap of Cen X-3 
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is "measured" as ~2-3 km. This is the first, albeit model dependent, measurement of the size of the 
accreting polar cap of an X-ray pulsar. 



5. Discussion of Low Frequency Behavior (< 100 Hz) 

A simple power law with a slope ~ —2 matches the broadband continuum from 3 — 100 Hz in the 
power density spectrum of Cen X-3 (see Figure [l], Figure ^| and Table 2). This general behavior has 
been previously reported for Cen X-3 and other HMXB sources (Pelloni and Hasinger 1990| ). 

The origin of this lower frequency noise may be the consequence of the physics controlling the loading 
of mass onto the magnetospheric field, perhaps modulated by intrinsic turbulence in the disk. In the 
case of Cen X-3, the mass flow exterior to the magnetosphere is almost certainly in an accretion disk. 
To fix ideas, suppose the boundary layer between the disk and the magnetosphere is Ray leigh- Taylor 
unstable, a mass entry mechanism previously treated in some detail for the case of spherical accretion 
by |Arons and Lea 1976| , |Arons and Lea 1980| , [Eisner and Lamb 1977| , [Eisner and Lamb 1984| . Some 



aspects of Rayleigh- Taylor instability in the disk/magnetosphere interaction have also been discussed 
by |Taam and Spruit 1990| . 

The thickness and height H = Rm{c ms /vK) of the disk's boundary layer set the size of an unstable 
mode (a "blob"). Then 5M ~ pH 3 is the mass of a blob, while in statistically steady accretion 
p = M /27rHR m VK (see |Arons 1993| for this result, in the relevant case when the plasma is diamagnetic, 



and |Chosh and Lamb 1979| when macroscopic dissipation is assumed to give the free access of the 



neutron star's magnetic field to the disk's plasma on microscopic scales.) Here R m ~ 4300 km is the 
magnetopause radius, evaluated assuming approximate corotation between the Keplerian disk's inner 
edge and the magnetosphere ( |Bildsten et al. 1997| ), vk ~ 5580 km/s is the Kepler velocity of the disk 



just outside the magnetopause, and c ms = (c 2 + v\Y^ 2 is the magnetosonic speed in the disk, with c s the 
sound speed and v a the Alfven speed determined by the disk's own magnetic field. Compton scattering 
of the hard X-rays from the star limits the sound speed to c s ~ 1000 km s -1 . As we shall show, in mass 
transfer by blobs, the disk's inner edge must have a strong magnetic field of its own, such that v\ 3> c 2 . 

The mass of a blob is 5M = pH 3 ~ MH 2 /2nR m VK, and the number of blobs contributing to the 
mass flow is « 2nR m /H. Therefore, in a statistically steady state, 

M^ WM ~^^Mf. (4, 
at riff t R m 

Here r is the number of magnetospheric free fall times t// = R m /vff(R m ) = (R 3 n /2GM*) 1 / 2 ~ 0.46 s 
required for a blob to form, decelerate and freely fall through a distance ~ H, thus detaching from the 
boundary layer. This detachment process is very ill understood. If one assumes the blob loses its orbital 
velocity instantaneously, then r = 4 (H/R m <C 1), while r ~ 3 if H/R m ~ 1. r is expected to be larger 
still, if one accounts for the mechanical and electromagnetic drag forces on the forming blob. 
Eliminating M from f| yields 

H c ms r 3 . . 

— — = —= —= ~ 2. (5) 

Rm v K v/2 y/2 



12 



Therefore, the inner edge of the disk must be reasonably thick, with iV^ofc ~ 2y / 27r/r ~ 3 blobs 
contributing at any one time to the mass flow, and the magnetic pressure in the disk must be large 
compared to the gas pressure. For Cen X-3's parameters, we find B 2 /47r ~ 230nkT gas . Note that we 
have neglected the inflationary effects of radiation pressure on the disk. Since Cen X-3's luminosity is 
~ 30 % of the Eddington luminosity, some of the disk's inflation might be supported radiatively rather 
than magnetically, but it is hard to avoid the conclusion that magnetic support makes a large and 
probably dominant contribution to the disk's inner structure. 

Such variability in the accretion flow will cause the size of the polar cap to vary. In the "kinematic" 
estimates of the size of the polar function of M (|Pringle and Rees 1972| , Pavidson and Ostrikeij 

1973n [Lamo", Pethick and Pines 19731) , 6 C oc {R*/R m ) 1/2 oc M 1/7 . Such a model might apply if a) 
macroscopic reconnection causes all dipole field lines which would have closed at radii r > R m to be 
open and b) blobs do not fall to radii substantially smaller than R m before their particles become frozen 
to magnetospheric field lines. If, on the other hand, blobs form diamagnetically, lose their angular 
momentum and fall freely before disrupting into particles which freeze to the magnetospheric field, then 
9 C oc M~ 2 / 7 (|A rons and Lea 1980|) . If blob formation causes the mass flux onto one pole to vary by as 
much as ^/l/(iV Wob /2), then in Cen X-3, the polar cap size might be varying by as much as 12 % (open 
field model) or 23 % (freely falling blob model) over the course of the three days of observation that 
went into forming the power spectrum. In principle, we should compare the average of several models 
similar to model A but with different values of 9 C to the high frequency fluctuation power, rather than 
the single "snapshots" shown in models A, B and C. However, the expected small range of variation in 
6 C suggests the snapshot comparison to be good enough, given the limitations in the simulations. 

Consideration of Kelvin Helmholtz mixing {e.g. |Burnard, Lea and Arons 1983| ); lead to similar 



conclusions. The KH effects most likely contribute to blob formation and deceleration, while Rayleigh- 
Taylor (RT) instabilities control the detachment of blobs from the disk's inner edge. For both effects, 
the basic mode size is H, and the mass transfer rate is controlled by gravitational free fall, so our 
order of magnitude estimates are invariant to the greater complexity of considering both RT and KH 
together. On the other hand, if the mass transfer were to be controlled solely by reconnection, with 
a mass loading velocity ~ ev^ with e perhaps on the order of 0.1 (see [Lazarian and Vishniac 1999 
for a recent investigation of rapid reconnection and full references on the subject). Repetition of our 
arguments lead to H/R m ~ 27r/e 1, which is impossible since the magnetic stresses would cause the 
disk to explode and shut off accretion. 

These estimates lead to several important conclusions. Only N b i ob ~ 3 "blobs" contribute to the 
mass flow onto each polar cap at one time, with the basic time scale for these slow fluctuations being 
rtff ~ 1.4 seconds. Thus, the observed fact that spectrum flattens at less than 3 Hz can be understood 
as the consequence of the fluctuations becoming large amplitude at these low frequencies. Possibly 
additional fluctuations are imposed on the structure of the blobs by intrinsic disk turbulence. Hawley 
has calculated models of disk turbulence which show a scale free power law dependence of the spatial 
transform of density, velocity and kinetic energy of the flow flHawIcy, Cammic and Balbus 1995| ). The 
necessary truncation of the power law on the largest scale corresponds to overturning velocity structure 
of order the inner disk scale height. These slowest components of the flow could naturally explain 
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the flattening of the temporal transform at low frequency (< 3 Hz). The power law form of the 
spatial transform would naturally imply a power law form of the temporal transform. The reality 
of these fluctuations depends on the assumption that the disk variability imposes a structure on the 
formation of blobs. Further, since the Rayleigh Taylor modes can include structure elongated along 
the circumference of the magnetosphere ("blob" azimuthal lengths between ~ H and ~ 2irR m ), lower 
frequency power might also be caused by modulations of the flow due to azimauthal structure in the 
instabilities, especially if one takes account of the time needed to get rid of the orbital velocity and 
allow a blob to go into free fall. Quantitative consideration of these more advanced issues is far beyond 
the scope of this paper. However, whatever mechanism modulates the accretion rate and accounts for 
the low frequency (3 — 100 Hz) power law with a slope ~ —2, its behavior is likely to extrapolate as a 
continuation of the same power law to higher frequency (> 100 Hz). Therefore the deviation from this 
power law at > 100 Hz that is observed in Cen X-3 is likely a different physical mechanism associated 
with the surface of the neutron star. 



6. Conclusions 

We have advanced the view that the fluctuations observed in the power spectrum at frequencies 
exceeding a few hundred Hz have their origin in the intrinsic photon bubble instability of the accretion 
flow at the stellar surface, rather than being the passive result of still higher frequency fluctuations 
in the mass loading into the magnet sophere. Emprically, this view is supported by the change in the 
power spectrum's structure above a few hundred Hz. From the point of view of polar cap plasma 
dynamics, the large amplitude fluctuations in the mass flux implied by the unstable mass entry model 
outlined in the previous section can be modeled as a steady flow, since the time variation scale of the 
mass flux (~ 0.5 sec for Cen X-3) is far longer than the intrinsic photon bubble instability time scale 
(milliseconds). Therefore, our radiation hydrodynamic calculations, which assume steady flow on the 
photon bubble time scale, are directly applicable. Our identification of the high frequency fluctuations 
in Cen X-3 (kHz QPO and broadband kHz continuum) with the photon bubbles (PBO and turbulence) 
is well supported by the fact that the four parameter radiation hydrodynamical simulations match the 
observations remarkably well in the applicable frequency range (> 100 Hz), with an adjustable polar 
cap size 6 C as the only free parameter of the fit. The mass of the neutron star, the magnetic field 
at its surface and the accretion rate are all constrained by other observations. Adjusting 9 C yields 
simulation calculations that qualitatively and approximately quantitatively explain eleven (of fifteen) 
free parameters of the best fit to the observed data. These include the form of the continuum power 
with a plateau, the high corner frequency and rolloff power law and the presence of two broad lines 
(QPO) at the appropriate frequencies, amplitudes and widths. The remaining four parameters of the 
fit to the observed data concern the low frequency (< 100 Hz) behavior of the source and the PCA 
deadtime corrections. The range of applicability of the calculations does not include behavior of the 
source for frequencies below ~ 100 Hz. The mass entry model outlined above suggests that these low 
frequency effects are related to modulation of the accretion rate due to the intrinsic instability of the 
mass entry mechanism, to which the X-ray emission from the surface passively responds. Thus, the 
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match of the models to the observations lends strong support that our discovery of kHz fluctuations in 
the X-rays emitted by Cen X-3 can be identified with our prior suggestion that PBO would be observed 
in this X-ray source ( [Klein et al. 1995a] ). 

The similarity between GRO J1744-28 and Cen X-3 is significant and strengthens the PBO identi- 
fication for GRO J1744-28 flKlein et al. 1996b|) . Both GRO J1744-28 and Cen X-3 show two QPO 
peaks and a strong continuum in their power density spectra ( [Zhang et al. 1996|) . The less complex 
pulse profile of GRO J 1744-28 (nearly a perfect sinusoid) suggests a larger polar cap than for Cen X-3 
in agreement with observations of lower frequency PBO/QPO. It is clear that the ~ 400 Hz separation 
in the kHz QPOs in Cen X-3 has nothing to do with rotation since this separation frequency is vastly 
different than the 0.21 Hz spin frequency. 

In contrast with the interpretation of GRO J 1744-28, the identification of photon bubble phenomena 
in CenX-3 is much firmer because of the naturally faster time scales extending to ~ 1800 Hz in the broad 
continuum. The observations reported here reaching to ~ 1800 Hz are the fastest known variations in 
the X-ray emission from an accreting compact object. We can not exclude the simultaneous presence of 
additional mechanisms for temporal variability above 100 Hz but PB models provide the only mechanism 
which yields quantitative estimates for both amplitudes and frequencies. 

Although the PB simulations include most of the accretion physics thought to be important near the 
surface of a neutron star undergoing optically thick accretion at the surface, there are some practical 
limitations. These models can currently only be integrated for a few tens of milliseconds, which limits a 
detailed comparison with averaged power density spectra over days of exposure in an observed source. 
Future X-ray missions with significantly larger area detectors and parallelization of our algorithms would 
circumvent these difficulties. 

Prior to the RXTE launch we predicted the observation of photon bubble phenomena including the 
specific existence of kHz fluctuations and some power in quasi-periodic oscillations identifiable as PBO 
in Cen X-3 ( [Klein et al. 1996a] ). We have now observed these fluctuations in Cen X-3 and have refined 
the model calculations to match the actual frequency dependence of the power density spectra including 
the presence of a significant continuum. A careful comparison of the predicted power density from the 
PB simulations with the observed power density spectrum constrains the size of the polar caps of Cen 
X-3 to a radius of ~2-3 km. This has permitted the first model-dependent measurement of the size of 
the accreting polar cap in an X-ray pulsar. 

The identification of PB phenomena in Cen X-3 opens the door to future X-ray observations with 
high resolution timing that will allow us to probe the physics of the surface of super-Eddington accreting 
X-ray pulsars. 
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8. Appendix A: PCA Deadtime Corrections 

In this appendix we describe the procedure for correcting the power density spectrum for deadtime 
effects in the PCA. These corrections are required and are important in those cases in which the Leahy 
normalized power density falls below the level of 2 that is expected if no deadtime effects are present. 
We use a frequency dependent model for an additive deadtime component of the power density of a 
form given by equation || This Leahy normalized deadtime model Pat(f) is suggested by equation (22) 
of |Zhang et al. 1995] and equation (2) of [Zhang et al. 19961 where f is the frequency. The functional 
form for the paralyzable and nonparalyzable deadtime models are nearly the same (compare Figures 1 
and 3 in Zhang et al. 1995 ). The practical distinction between the two models is insignificant below 
10 kHz. The effect of discrete time binning is also unimportant below 10 kHz for the data shown in this 
appendix. The value of td is about 10/is and is a fixed property of the function of each of the five modules 
that comprise the PCA. The value of t v i e is a selectable parameter of the flight configuration of the PCA 
and is ~ 61/xs (level 1) or ~ 142/xs (level 2) for all data shown in this appendix. These selectable values 
of t v i e are averages of separate ground calibration measurements for the five PCUs of the PCA. The 
variations in these values of t v \ e are of order 1% which is more that sufficient for the purpose of evaluating 
deadtime corrections. The variables r G and r in are the actual counting rate and the incident counting 
rate (assuming no deadtime correction) respectively. The variable r v \ e is the total rate of VLE (very 
large event) due to background particles and pileup of X-rays. These rates correspond to the values for 
one of the five PCA modules. Clearly the values of r a and r in are complex time dependent functions as 
revealed by the actual observations of r a . Equation |6] is suggested by the simple situation in which both 
r Q and rj n are constant. With further simplifying assumptions that r Q rs Ti n and r Q <C t~[ , we can 
express the values of A and B as approximate functions of the rates r Q and r v i e (see equations |7] and §). 
We will show that the overall deadtime correction Pdt{f) can be successfully represented by fixed values 
of A and B for actual observations. Note that our proposed procedure for estimating the deadtime 
corrections is to vary the values of A and B as additional free parameters of any model that one might 
use to fit a power density spectrum. The approximate relationship of A and B to any rates (r , Tj n or 
r v i e ) as given by equations [7] and [8] serves only as a verification for sensible values of A and B since these 
equations are based on a model that assumes that the rates are constant when in reality these rates 
are highly variable. Parameter C and third term of equation^ is the additive component of background 
present in all PCA observations. Its functional form was determined by approximately (within about 
1%) fitting the empirical power density spectrum of a deep exposure to background with no signicant 
sources present QMorgan 1999 ). This frequency dependent function shows clear minima at 2000, 6000 
and 10000 Hz and maxima at 0, 4000, and 8000 Hz. These features are likely produced by background 
events in channels 0-7. The total count rate from the background for GX 17+2, Cyg X-l and Cen X-3 
in channels 0-17 is ~ 30 s _1 . Since the background rate is not separately measured in the presence of 
a bright source there is no firm estimate of the background rate during an typical observation. The 
conservative procedure is to include an additive term in the fit to any observe power spectra with an 
amplitude proportional to the parameter C. This background term is not really a deadtime correction 
but is encluded in the procedure since it is a model for the frequency dependent form of a needed 
correction. In the worst case an error in the determination of the values of A, B and C would correlate 
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with some real component of the power density spectrum of an observed source. This possibility always 
leads to an underestimate of the actual power density spectrum of the source. Therefore this procedure 
which freely varies the parameters A, B and C is a conservative approach for estimating any components 
of the power density spectrum which are not real components of the observed target. 
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In Figure [5] we show the power density spectrum of GX 17+2 over a frequency range that extends to 
~ 50 kHz (|M organ 1999| ). The curve shows the Leahy normalized power density on a linear scale near 
the expected value of 2 as a function of the logarithm of the frequency. The value of the power density 
falls below 2 for frequencies greater than ~ 100 Hz, clearly indicating that large deadtime corrections 
are necessary to determine the intrinsic component of the power density. The solid curve is the deadtime 
component of the power density spectrum which approximately matches above ~ 2000 Hz as defined by 
equation |] with values of A and B equal to 0.0597 and 0.0101 respectively. The value of C is consistent 
with indicating the low level of any background effects for source as bright as GX 17+2. The match 
of the deadtime model to the observed data is excellent above ~ 2000 Hz as indicated by the one sigma 
error bars (dotted stairstep curves). In particular the broad dip in the power density near ~ 5000 Hz is 
fully explained as a deadtime correction due to VLE effects. The lower dashed curve shows the deadtime 
model if B = with no VLE correction. The upper dashed curve shows the deadtime model if the 
VLE correction is arbitrarily increased by a factor of six. This extreme form of the deadtime model 
shows a feature at ~ 9 kHz which is clearly in excess of the power density level near that frequency. The 
correct deadtime model (solid curve) matches the slight "bump" in the power density near ~ 9 kHz. 
Figure |6| is an expanded version of the data shown in figure |^. The purpose of this figure is to clearly 
show the match of the deadtime model to the slight bump in the power density near ~ 9 kHz. This 
feature corresponds to the second peak of the 



sin(x) 



function for the VLE correction to the deadtime 
model and therefore depends on the ground calibrated value of t v \ e . The best values of A and B are 
approximately consistent with values of r Q ~ 1500s -1 and r v \ e ~ 120 s -1 for one of the five PCA 
modules. The implied total rates in the PCA are approximately consistent with the actual measured 
rates. The most conservative procedure for using this deadtime model is to allow the values of A and 
B to vary in order to match the high frequency component of the power density spectrum which is 
likely to be dominated by deadtime effects. One can then derive approximate values of r a and r v[e by 
inverting equations [7| and § and then checking that these values are a rough match to the typical 
measured rates. With this interpretation of the deadtime correction Figure [5] shows the presence of 
the kHz QPO at ~ 600 Hz previously discovered by |Wijnands et at. 1997] but no significant intrinsic 
power above ~ 1000 Hz. This deadtime corrected interpretation of the power density spectrum of GX 



18 



17+2 shows for the first time that there is very little continuum power on either side of the kHz QPO 
feature. In summary, this power density spectrum of GX 17+2 indicates that a simple deadtime model 
with two free parameters A and B successfully explains the frequency dependent form of the deadtime 
correction to an accuracy of ~ 1% of the depression of the Leahy normalization 2 which is expected for 
pure Poisson noise with no deadtime correction. 

Since our goal is to demonstrate that Cen X-3 shows excess power density in the range from 300 — 
2000 Hz above the level of the deadtime correction, we seek an independent source of sufficient counting 
rate other that GX 17+2 which has no known intrinsic power in this frequency band. Figure [7] shows 
the power density spectrum of Cyg X-l over a frequency range that extends to ~ 6000 Hz. The data for 
the curve shown in Figure [7| were provided by |Gui 1999] and show the Leahy normalized power density 
on a linear scale near the expected value of 2 as a function of the logarithm of the frequency. The value 
of the power density falls below 2 for frequencies greater than ~ 100 Hz clearly indicating that deadtime 
corrections are necessary to determine the intrinsic component of the power density. The solid curve is 
an approximate fit of the deadtime corrected power density above ~ 300 Hz as given by equation || The 
approximate best values of A and B are 0.05 and 0.01 respectively which corresponds to mean counting 
rate r ~ 1250 and r v i e ~ 140 for one of the five PCA modules. The implied total rates in the 
PCA are approximately consistent with the actual measured rates and are similar to the rates for GX 
17+2. Again the background term (parameter C) is not needed because of the brightness of CYG X-l. 
The match of the deadtime model to the observed data is excellent above ~ 300 Hz. The lower dashed 
curve shows the deadtime model if B = with no VLE correction. The upper dashed curve shows 
the deadtime model if the VLE correction is arbitrarily increased by a factor of three which is clearly 
inconsistent with the power density spectrum. The deadtime corrected interpretation of this spectrum 
indicates the significant detection of intrinsic power from Cyg X-l up to ~ 280 Hz. This interpretation 
extends and is consistent with prior published observations of Cyg X-l ( Pui et al. 1997] and |Nowak e\ 
al. 19|99|). For frequencies greater than ~ 300 Hz the spectrum is consistent with the deadtime model 
including the slight drop in power near ~ 2000 Hz. Since Cyg X-l is a few times brighter than Cen 
X-3, these deadtime corrections will be greater than those for Cen X-3. The results for the Cyg X-l 
data show that deadtime effects from 300 — 4000 Hz are smooth and correctable with no unexplained 
features. These data from Cyg X-l were obtained in the similar PCA mode (SB_125us_0_17_ls)) as the 
data from Cen X-3 reported in this paper (SB_62us_0_17_ls). In particular the energy range for the two 
data sets is identical. One could successfully use a scaled version of the spectrum of Cyg X-l between 
300 — 4000 Hz as the deadtime model for the Cen X-3 observations. 

The stair step solid curve in Figure || shows the power density of data from Cen X-3 over a frequency 
range ~ 50 Hz to ~ 10 kHz. The stair step dotted curves above and below are the one sigma errors of 
the measured power density. These data are the same as shown in both figures [l] and [| and are plotted 
here to clearly show the level of the necessary deadtime correction. 

The curve shows the Leahy normalized power density on a linear scale near the expected value of 2 as 
a function of the logarithm of the frequency. The value of the power density falls below 2 for frequencies 
greater than ~ 50 Hz clearly indicating that deadtime corrections are necessary to determine the intrinsic 
component of the power density. The upper solid curve is the best fit model (15 parameters) to the total 
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power density as a function of frequency. The range and best fit values of each of the 15 parameters is 
shown in Table 2. Three of the parameters are the components A, B and C of the deadtime model in 
the form given by equations and |S|. The lower solid curve shows the deadtime model for the best 
fit values of A, B and C (see Table 2). These data for Cen X-3 were obtained with a setting for t v \ e 
of ~ 61/is (level 1) in contrast a higher setting of 142/is(level 2) for the GX 17+2 and Cyg X-l data. 
The lower setting decreases any VLE deadtime effect and was set to this lower value for the Cen X-3 
observations precisely to decrease any VLE effect. The values of A and B are used to estimate typical 
values for r D (~ 310 and r v i e (~ 40 which roughly agree with the observed average rates. 
The match of the deadtime model to the observed data is excellent above ~ 2 kHz. The lower dashed 
curve shows the deadtime model if B = with no VLE correction. The upper dashed curve shows the 
deadtime model if the VLE correction is arbitrarily increased by a factor of six. Such a model roughly 
matches the level of the power density in the frequency range from 100 — 1000 Hz but is completely 
inconsistent with the power density from 1 — 10 kHz. Clearly the excess power density in the frequency 
range 100 — 2000 Hz can not be explained as a VLE deadtime effect. The dotted continuous curves just 
above and below the lower solid curve show the affect on the deadtime model if the VLE correction 
factor (parameter B) is varied from a low to a high value corresponding to the 50% confidence range. 
A correct VLE deadtime correction does effect the determination of the intrinsic level of the continuum 
power in the 100 — 1000 Hz range. However the continuum component (see dashed-dotted curve) is well 
above the level of the deadtime component (lower solid curve). The lower short- dash- log- dash curve is 
the upper limit to the background as determine by the free parameter C. It shows a drop at 2000 Hz 
followed by a small peak at 4000 Hz which reachs but does not exceed the measured power density at 
4000 Hz. This upper limit for C greatly exceeds any expected level of the background which is likely 
at least full order of magitude lower. Even so, a value of C that is increased by a factor 8 (see the 
upper short-dash-long-dash curve) does roughly match the overall shape of the power density spectrum 
from 100-2000 Hz. This artificially amplified value of C would create a strong peak at 4000 Hz which 
is well above the measured level of the power density near 4000 Hz. Estimates for the contribution of 
the background are also shown in Figure |l|. An additive background contribution is not a significant 
component of the power density spectrum of Cen X-3. 
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TABLE 1 
Parameters of Cen X-3 



Lower Value 



Upper Value 



Mass of NS 


M„s 


1.4 


1.3 


1.5 


solar mass 


Radius of NS 


Una 


10.0 


9.0 


17.0 


km 


Cyclotron Line 


E C 


28.5 


27.5 


29.5 


kcV 


Magnetic Field 


B 


3.2 


3.1 


3.3 


10 12 Gauss 


Distance 


D 


10.0 


9.0 


11.0 


kpc 


X-ray Luminosity 


L, r 


9.4 


9.0 


9.8 


10 37 ergs _1 


Polar Cap Radius 


6 C 


0.25 


0.1 


0.4 


radians 



PnthirV 



iantanffoln nf nl 1 QQfj 



at 



rlutchings et al. 19791 



(this paper) 
(this paper) 
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-2 





Log Frequency(Hz) 



Fig. 1. — The curve is the observed average power density in units log x 



Hz 



-i 



as a function 



510 

of log 10 [frequency (Hz)]. Note the presence of the first and second harmonics of the pulse (0.21 and 
0.42 Hz) due to the rotation of the neutron star, the power law continuum below ~ 100 Hz with a 
slope of ~ —2 (dashed line), the flattening of the curve near 100 Hz, the two kHz QPO peaks, and the 
continuum power between 100 — 2500 Hz. This intrinsic power density spectrum has been corrected 
for the effects of Poisson noise including the slight frequency dependence caused by deadtime effects 
of the PCA (see appendix A). The dotted curve shows an estimate of the upper limit to the level of 
the deadtime corrected Poisson level. Power density at high frequency that falls below the dotted line 
is not confirmed as intrinsic to Cen X-3. The long-dashed curve just below the dotted curve shows 
the upper limit to any additive contribution from the background. This non-source component of the 
power density spectrum can not be much larger that shown here otherwise the peak at 4000 Hz would 
be detectable. The lower dot-dashed curve is the most likely level of the background based on the 
correction procedure (see Appendix A). The actual level of the background is likely another order of 
magitude smaller than the dot-dashed curve. 



23 



-2 



2 
Log Frequency(Hz) 



Fig. 2. — The curve of measured data (solid stair step) is the same as in Figure |l], but has now been 
multiplied by frequency. Dotted lines are plotted above and below the data at la error levels. Dashed 
lines are plotted horizontally at levels of 2, 5 and 100 % C^^j- Note that the pulsed fraction is at a 
level of several tens of %. The two kHz QPO peaks and the continuum power above 100 Hz have a total 
power bounded between 2 and 5% rms. The best fit to the data plotted as a continuous solid curve 
tracks the measurements. 
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8 



3 4 
Log Frequency(Hz) 



8 



3 

Log Frequenc 



Fig. 3. — The left panel is a log-log plot of power density spectra in units The right panel is a smoothed 
representation of the data in the left panel and is more appropriate for direct comparison with the 
observed data. The vertical scale has been shifted up by 2 and 4 decades for models B and C respectively 
for illustrative purposes. Note that all three models show evidence of PBOs and PB continuum power. 
Lines that indicate power laws of the form f~ 2 (dotted) and (dashed) are indicated in the right 

panel. These slopes at high frequency are roughly consistent with the results that we reported previously 
i 
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a. 



o —7 
J 1 



Log Frequency(Hz) 



Fig. 4. — The upper curve in each panel is the observed power density of Cen X-3 from Figure |2|. 
The error bars are shown to facilitate a quantitative comparison to the theoretical models. Beyond 
~ 1800 Hz the 2a upper limit is shown as a dotted curve. Also the low frequency power law fitted to 
the observed data has been subtracted since this component of the observed spectrum is not included 
in the physics incorporated into the theoretical models. The power density spectra of the three models 
A, B and C are shown on the same absolute scale. These are the smoothed versions of the spectra 
appropriate for direct comparison to the observed spectrum of Cen X-3. The three models have been 
correctly red shifted by the appropriate amount for the surface of a neutron star. 
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1.92 



3 

Log Frequency 



Fig. 5. — The power density of data from GX 17+2 over a frequency range that extends to ~ 50 kHz. 
The curve shows the Leahy normalized power density on a linear scale near the expected value of 2 as a 
function of the logarithm of the frequency. The value of the power density falls below 2 for frequencies 
greater than ~ 100 Hz clearly indicating that large deadtime corrections are necessary to determine 
the intrinsic component of the power density. The solid curve is an approximate fit of the deadtime 
corrected power density above ~ 2000 Hz as given by equations |6], [7| and §. The match of the deadtime 
model to the observed data is excellent above ~ 2000 Hz. The lower dashed curve shows the deadtime 
model if B = with no VLE correction. The upper dashed curve shows the deadtime model if the 
VLE correction is arbitrarily increased by a factor of six. This extreme form of the deadtime model 
shows a feature at ~ 9 kHz which is clearly in excess of the power density level of the data. The correct 
deadtime model (solid curve) matches the slight "bump" in the power density near ~ 9 kHz. 
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Log Frequency 



4 



Fig. 6. — This figure is an expanded version of the data shown in figure || The purpose of this figure is 
to clearly show the match of the deadtime model to the slight bump in the power density at ~ 9 kHz. 



This feature corresponds to the second peak of the 
deadtime model. 



sin(x) 
x 



function for the VLE correction to the 



TABLE 2 
Fitted Model Parameters 



Parameter 


Variable 


Value 


Lower Value 


Upper Value 


Units 


Low Corner Frequency 


Fi 


63 


46 


76 


Hz 


High Corner Frequency 


F h 


1228 


923 


1662 


Hz 


PD(F,) 


Pi 


2.50 


2.05 


2.95 


10- 6 (-£2g_) 2 Hz" 1 


Low PD Slope 


Si 


-2.15 


-3.89 


-1.04 


dimcnsionlcss 


Mid PD Slope 


s m 


-0.27 


-0.37 


-0.18 


dimcnsionlcss 


High PD Slope 


s h 


-2.61 


-8.1 


-1.5 


dimcnsionlcss 


Frequency of QPO 1 


Fi 


331 


260 


407 


Hz 


Amplitude of QPO 1 


Ai 


3.1 


2.1 


3.9 


% rms 


Equiv. width of QPO 1 


Ai 


712 


329 


1106 


Hz 


Width of QPO 1 (1 a) 


Wi 


89 


48 


213 


Hz 


Width of QPO 1 (FWHM) 


Wi 


214 


115 


511 


Hz 


Q of QPO 1 (FWHM) 


Fi/Wi 


1.5 


1.2 


1.9 


dimcnsionlcss 


Frequency of QPO 2 


F 2 


761 


671 


849 


Hz 


Amplitude of QPO 2 


A 2 


3.1 


2.0 


4.0 


% rms 


Equiv. width of QPO 2 


A 2 


969 


426 


1534 


Hz 


Width of QPO 2 (1 a) 


W 2 


103 


36 


262 


Hz 


Width of QPO 2 (FWHM) 


w 2 


247 


86 


629 


Hz 


Q of QPO 2 (FWHM) 


F 2 /W 2 


3.1 


2.7 


3.4 


dimcnsionlcss 


Deadtime Factor 


A 


0.01304 


0.01285 


0.01328 


dimcnsionlcss 


VLE Deadtime Factor- 


B 


0.00071 


0.00048 


0.00095 


dimcnsionlcss 


Background 


c 


0.00011 


-0.00032 


0.00055 


dimcnsionlcss 
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1.95 



Log Frequency 



Fig. 7. — The power density of data from Cyg X-l over a frequency range that extends to ~ 6000 Hz. 
The curve shows the Leahy normalized power density on a linear scale near the expected value of 
2 at as a function of the logarithm of the frequency. The value of the power density falls below 2 
for frequencies greater than ~ 100 Hz clearly indicating that deadtime corrections are necessary to 
determine the intrinsic component of the power density. The solid curve is an approximate fit of the 
deadtime corrected power density above ~ 300 Hz as given by equations [| |7] and |[ The match of the 
deadtime model to the observed data is excellent above ~ 300 Hz. The lower dashed curve shows the 
deadtime model if B = with no VLE correction. The upper dashed curve shows the deadtime model 
if the VLE correction is arbitrarily increased by a factor of three. 
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1.985 

1 



2 3 
Log Frequency(Hz) 



Fig. 8. — The power density of data from Cen X-3 over a frequency range ~ 50 Hz to ~ 10 kHz 
is the stair step solid curve. The stair step dotted curves above and below are the one sigma errors 
of the measured power density. These data are the same as shown in both figures [I] and |2| and are 
plotted to clearly show the level of the necessary deadtime correction. The curve shows the Leahy 
normalized power density on a linear scale near the expected value of 2 as a function of the logarithm 
of the frequency. The upper solid curve is the best fit model (15 parameters) to the total power density 
as a function of frequency. The lower solid curve is the best fit of the deadtime correction to the power 
density. The match of the deadtime model to the observed data is excellent above ~ 2 kHz. The lower 
dashed curve shows the deadtime model if B = with no VLE correction. The upper dashed curve 
shows the deadtime model if the VLE correction is arbitrarily increased by a factor of six. The dotted 
continuous curves just above and below the lower solid curve show the effect on the deadtime model if 
the VLE correction factor (parameter B) is varied from a low to a high value that corresponds to a 50% 
confidence range. The continuum component (see dashed-dotted curve) is well above the level of the 
deadtime component (lower solid curve). 
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1 



2 

Log Frequency(Hz) 



Fig. 9. — The power density residuals of data from Cen X-3 over a frequency range ~ 50 Hz to ~ 10 kHz 
is the stair step solid curve. The model minus data residuals are plotted in units of the local measured 
one sigma error for each data point. These model and data are the same as shown in figure ||. The 
model approximately uniformly fits the data above ~ 100 Hz and therefore all parameters related to 
the kHz QPO and continuum are formally well defined. The power law model below ~ 100 Hz is clearly 
not sufficiently complex to completely model the behavior of Cen X-3 at low frequency. This difficulty 
with the model is amplified by the excellent signal to noise below ~ 100 Hz. 
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